function K = covQPeriodic(hyp, x, z, i)

% Stationary covariance function for a smooth periodic function, with period p:
%
% k(x,y) = sf2 * exp( -2*sin^2(pi*||x-y||/p )/ell^2 )
%
% where the hyperparameters are:
%
% hyp = [ log(ell)
%         log(p)
%         log(sqrt(sf2)) ]
%
% Copyright (c) by Marco A F Pimentel, same format as Carl Edward Rasmussen 
% and Hannes Nickisch, 2012-10-22.
%
% See also COVFUNCTIONS.M.

if nargin<2, K = '5'; return; end                  % report number of parameters
if nargin<3, z = []; end                                   % make sure, z exists
xeqz = numel(z)==0; dg = strcmp(z,'diag') && numel(z)>0;        % determine mode

% n = size(x,1);
ell = exp(hyp(1));
p   = exp(hyp(2));
sf2 = exp(2*hyp(3));
pm  = exp(hyp(4));
m   = exp(hyp(5));

% precompute distances
if dg                                                               % vector kxx
  K1 = zeros(size(x,1),1);
  K2 = zeros(size(x,1),1);
else
  if xeqz                                                 % symmetric matrix Kxx
    K1 = sqrt(sq_dist(x'));
    K2 = sqrt(sq_dist(sin(2*pi*x'/pm)));
    K3 = sqrt(sq_dist(2*pi/pm*x'.*cos(2*pi*x'/pm)));
  else                                                   % cross covariances Kxz
    K1 = sqrt(sq_dist(x',z'));
    K2 = sqrt(sq_dist(cos(2*pi*x'/pm),cos(2*pi*z'/pm)));
    K3 = sqrt(sq_dist(2*pi/pm*x'.*cos(2*pi*x'/pm),2*pi/pm*z'.*cos(2*pi*z'/pm)));
  end
end

K1  = pi*K1/p;
K2  = m*K2/2;
if nargin<4                                                        % covariances
    K = sin(K1+K2)/ell; K = K.*K; K =   sf2*exp(-2*K);
else                                                               % derivatives
  if i==1
    K = sin(K1+K2)/ell; K = K.*K; K = 4*sf2*exp(-2*K).*K;
  elseif i==2
    R = sin(K1+K2)/ell; K = 4*sf2/ell*exp(-2*R.*R).*R.*cos(K1+K2).*K1;
  elseif i==3
    K = sin(K1+K2)/ell; K = K.*K; K = 2*sf2*exp(-2*K);
  elseif i==4
    R = sin(K1+K2)/ell; K = 2*sf2/ell*exp(-2*R.*R).*R.*cos(K1+K2).*K3;
  elseif i==5
    R = sin(K1+K2)/ell; K = -4*sf2/ell*exp(-2*R.*R).*R.*cos(K1+K2).*K2;
  else
    error('Unknown hyperparameter')
  end
end